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In recent work a deterministic and time-reversible boundary thermostat called thermostating 
' by deterministic scattering has been introduced for the periodic Lorentz gas [Phys. Rev. Lett. 84, 

4268 (2000)]. Here we assess the nonlinear properties of this new dynamical system by numerically 
(N . calculating its Lyapunov exponents. Based on a revised method for computing Lyapunov exponents, 

^ l' which employs periodic orthonormalization with a constraint, we present results for the Lyapunov 

O . exponents and related quantities in equilibrium and nonequilibrium. Finally, we check whether 

^0 ' we obtain the same relations between quantities characterizing the microscopic chaotic dynamics 

and quantities characterizing macroscopic transport as obtained for conventional deterministic and 
time-reversible bulk thermostats. 
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I. INTRODUCTION 



The investigation of nonequilibrium transport processes in many-particle systems generally requires to model the 
y—^ ' interaction between a particle and a thermal reservoir. A common approach for such a modeling are deterministic and 
I> ' time-reversible thermostats Conventional types of them, such as the Gaussian and the Nose-Hoover thermostat, 

are based on introducing a momentum dependent friction coefBcient into the microscopic equations of motion [^|-^. 
Though the microscopic equations of motion of these systems are time-reversible the macroscopic dynamics is irre- 
versible in nonequilibrium leading to momentum and energy fluxes with well-defined transport coefficients ||^,^,^ Q , 
which appears to be a paradox. However, investigations of the microscopic dynamics with methods from dynamical 
system theory could resolve this paradox by showing that the microscopic dynamics is nonlinear and highly unstable 
' and leads to a phase space volume contraction onto a fractal attractor 1 1 1 ijjl^ . From the analysis of con- 

, ventional thermostats, further relations between quantities characterizing the microscopic dynamics and quantities 
. 1-H ■ characterizing macroscopic transport could be established. At the heart of such relations there is an identity between 
(— I I phase space volume contraction and thermodynamic entropy production. On the basis of this identity the Lyapunov 
• • . exponents could be related to the transport coefficients of a system, which has been formulated as the Lyapunov sum 
.>: rule 001^]. 

. These characteristic features of thermostated many-particle systems have been recovered for specific one-particle 
■ systems, the Gaussian thermostated periodic Lorentz gas fll| , |l5| , |l8| , |l9| , ^ll -p^ and the Nose-Hoover thermostated 
5^ ', periodic Lorentz gas . The periodic Lorentz gas consists of a particle that moves through a triangular lattice of hard 
disks and is elastically reflected at each disk collision. It serves as a standard model in the field of chaos and transport 
@|6|. The advantage of a one-particle system is that it reflects more strongly and transparently the nonequilibrium 
properties induced by a thermostat. For this reason the Lorentz gas appears as an appropriate tool to compare the 
properties of nonequilibrium steady states obtained from different deterministic and time-reversible thermostating 
mechanisms. The study of different models describing the interaction between particles and thermal reservoir and 
the identification of their common properties is crucial to obtain a general characterization of nonequilibrium steady 
states. 
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To investigate whether the nonequihbrium properties of conventional deterministic and time-reversible thermostats 
are of general validity, or just characterize these specific types of systems, an alternative deterministic and time- 
reversible thermostat called thermostating by deterministic scattering has been introduced for the periodic Lorentz 
gas [§,0. This thermostat is based on specifically modeling the energy transfer related to a microscopic collision 
process between particle and disk, where the disk mimics a thermal reservoir with infinitely many degrees of freedom. 
In nonequilibrium under an external electric field this mechanism leads to an on average constant kinetic energy of the 
particle resulting in a nonequilibrium steady state. Furthermore, the phase space volume contracts onto an attractor 
similar to the multifractal attractor found for the Gaussian thermostated Lorentz gas. However, differences appear in 
the bifurcation diagram and in the field dependence of the conductivity. This alternative thermostat has later been 
applied to a heat and shear flow . 

In this work we focus on the microscopic properties of thermostating by deterministic scattering in the periodic 
Lorentz gas by numerically calculating the Lyapunov exponents. As quantities from dynamical systems theory, 
Lyapunov exponents allow a detailed characterization of the microscopic stability. In particular, they will enable us 
to check the general validity of relations between quantities from dynamical systems theory and statistical mechanics 
as obtained for conventional deterministic and time-reversible thermostats. We first explain the algorithm to calculate 
the Lyapunov exponents for the Lorentz gas as thermostated by deterministic scattering. Numerical computations 
will then show that the standard Gram-Schmidt orthonormalization has to be modified resulting in a variant of this 
method called constraint orthonormalization. Beside the results for the Lyapunov exponents we present results for 
the Kaplan- Yorke dimension and for the phase space volume contraction. We compare these results as obtained for 
our model with the results as known for the Gaussian thermostated Lorentz gas [^^, and with results for a heat 
and shear flow thermostated by deterministic scattering |^^. Finally, we check whether the phase space volume 
contraction is equal to the thermodynamic entropy production and whether the Lyapunov sum rule holds for our 
mechanism. 



II. ALGORITHM FOR THE CALCULATION OF THE LYAPUNOV EXPONENTS 



In a smooth d-dimcnsional system the equations of motion for a phase space vector F, 

r = F(r) , 

and the corresponding equations of motion for d tangent vectors 5T — (Sr, d-v), 



dF 



are integrated to obtain d Lyapunov exponents 



The Lyapunov exponents are a measure to characterize the stability of the dynamics | |3l| , |3^ . The maximal Lyapunov 
exponent A measures the maximal exponential divergence of two initially neighboring points ^r(O). However, during 
the time evolution every tangent vector will move into the fastest growing direction due to the instability of the 
dynamics. All these vectors will thus become indistinguishable and their norm will diverge. The algorithm of 
Benettin avoids this problem by a periodic Gram-Schmidt reorthonormalization of the tangent vectors thus enabling 
to compute the full spectrum of Lyapunov exponents associated to the d-dimensional phase space ||3^ , |3^ . 

In the periodic Lorentz gas the time-continuous flow describing the dynamics of a phase space volume vector T in 
the bulk is interrupted by a time-discrete map M describing the transformation of T at the moment of a collision, 

r' = M(r) . 

Dellago and coworkers have developed an algorithm to calculate the Lyapunov exponents for particle systems with 
hard sphere interactions and applied it to the Gaussian thermostated Lorentz gas ^ , p9|j35| . Here the tangent vectors 
are transformed at the moment of a collision according to the following rule p3| : 



dM 

^F(r)-F(M(r)) 



(2) 
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Eq. (g) is valid for arbitrary systems composed of a flow F and a time-discrete map M. It takes into account that a 
trajectory and a satellite trajectory collide with the disk at different space points and with a time delay 5tc, 



5Tr — 



(v,n) 



where n is the unit vector perpendicular to the surface at the collision point. 

Before we establish the equations of motions for the tangent vectors of the Lorentz gas as thermostated by deter- 
ministic scattering we briefly summarize the full equations of motion of Refs. UJ2^ for a particle described by the 
phase space vector F = (r,v). In the bulk F evolves according to 



r 

V = £< 



Vo 



(3) 



where e is an external electric field of strength e ~ \e\ generating a nonequilibrium situation. The basic idea of 
thermostating by deterministic scattering is now that at a collision energy is transfered such that the resulting velocity 
distribution for the particle is canonical in e quil ibrium. In a way, it results in a deterministic and time-reversible 
formulation of stochastic boundary conditions p?! , ^ . For this purpose the collision rules have been defined as follows: 
The velocity of the particle and its direction of flight are changed at a collision with the disk according to 



where 7 is the angle of incidence, X{'j) — sin I7I, B is the baker map and 



Yooiv) 



2 -i.^ V 
—ve~ + erf(— =) 



(4) 



(5) 



with T as a parameter corresponding to the temperature of the particle T =< > /2 at e = in equilibrium. The 
geometry of the periodic Lorentz gas and the relevant variables are shown in Fig. To ensure that the system is 
time-reversible, the forward baker map B acts if < 7 < 7r/2, and B is replaced by its inverse B"-*- if — 7r/2 < 7 < 0. 
To avoid any symmetry breaking induced by this combination of forward and backward baker map, we alternate their 
application in 7 with respect to the position /3 of the colliding particle on the circumference. For the spacing between 
two neighboring disks with the radius i? = 1 we choose, following the literature 21|, vj ~ 0.2361. Investigating 
this system in nonequilibrium by switching on an external electric field £ > leads to a nonequilibrium steady state 
with on average constant kinetic energy of the particle < >= const., i. e. the system is thermostated. 
The equations of motion for the tangent vectors ST in the bulk can now be derived from Eqs. as 



Sr' 



1 t 
1 



dr 



The transformation rules for the tangent vectors at the moment of a collision are obtained by inserting the collision 
rules for the phase space vector T of Eq. (^) into Eq. (|^), 



where 



1 
A B 



5r 
Sv 



1 
A B 



(6) 



dr-i 



and 



The components of the submatrix A read 

i-hi + l] 
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dvy' 
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dry 




dvy' 




dry 





[hi 



1] 



and the components of the submatrix B are 
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dvx 



y1 



I 



ftl + VxVy /l2 



dVx 

dVy 

9Vy' 



'^x ^y , , 1 1 

= ^/ll + VyVx «2 

= + VyVy /l2 



with 



hi = mi 



V (r, v') 



h2 = m2-pre' 



{mi, 7712} arc the slopes of the baker map {2, 0.5}, or the slopes of the inverse baker map {0.5, 2}, respectively. 

In the periodic Lorentz gas as thermostated by deterministic scattering the dynamics of four orthonormal tangent 
vectors has to be investigated to obtain four Lyapunov exponents which completely characterize the stability in the 
four-dimensional phase space. Before we present our results for the Lyapunov spectrum we wish to derive explicit 
expressions for two other interesting quantities. 



III. PHASE SPACE VOLUME CONTRACTION 



The phase space volume contraction P is equal to the sum of the Lyapunov exponents |32||3g]. 



(7) 



In the periodic Lorentz gas as thermostated by deterministic scattering only the change of a phase space volume 
element 5T at the moment of a collision as described by Eq. (||) contributes to P. The mean exponential rate of the 
phase space volume contraction P can then be calculated according to 



where 



d5T' 




dSV 





dSr' 







ddr ddv 

The partial derivatives of Eq. (||) read 



dSr 



d5rx 
dSry' 
dSrx 



1 + (Vx' - Vx) 



Vx) 



(v,n) 



(v,n) 



P =< In 



dST' 



dsr 



/ dSrx' dSvy' dSvx' d5ry'\^ ( dSvx dSvy' dSvx' dSvy' 
\ dSvx dSvy dSvy dSvx 



dSvx 06 Vy 



dSrx' 
dSvy 
dSr 



dSvy dSvx 



ddr 



= 1 + (Vy' - Vy) 



(v,n) 



(8) 



dSvx' 
dSvx 
dSvy' 
d6vx 



00 



B 



10 



Inserting these expressions into Eq. (|8|) yields 

P = 



< v'^ > - < > 

2r 



dSvx 
dSvy 
dSvy' 
dSvy 



= B, 



01 



Bi 



1 • 



(9) 



The phase space volume contraction depends thus only on the average transfer of kinetic energy to the reservoir. 
Equation is valid in equilibrium as well as in nonequilibrium. An analogous result has been obtained for collisions 
with a flat wall in a heat and shear flow thermostated by deterministic scattering p^ ], however, in case of shear the 
expression for P turned out to be more complicated. 
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If P < the phase space volume typically contracts onto a fractal attractor in the driven periodic Lorentz gas 
p] , P^Jl6| . The geometric properties of the attractor can be related to the Lyapunov exponents by the Kaplan- Yorke 
conjecture, Dky = Di. Here Di is the information dimension |^ and Dky is the Kaplan- Yorke dimension defined 
by 

Dky = 3 + > (10) 

where the Xi are ordered by magnitude, Ai > A2 > . . ., and j is the largest integer for which X]i=i 0. 

IV. THERMODYNAMIC ENTROPY PRODUCTION AND RESERVOIR TEMPERATURE 

The macroscopic properties of nonequilibriuni steady states can be characterized by quantities from thermodynamics 
and statistical physics. In this work we want to check whether we can relate the thermodynamic entropy production 
dS, 

dS='^, (11) 

to the phase space volume contraction. To calculate the thermodynamic entropy production for thermostating by 
deterministic scattering we have to calculate the temperature of the reservoir Tr in nonequilibrium. As discussed in 
Ref. p7[ |, in nonequilibrium the temperature related to the particle, or respectively the temperature in the bulk Tb 
defined via equipartitioning of energy, is greater than the parametric temperature T in Eq. and increases with 
the field strength. Moreover, Tb is inhomogeneously distributed in the bulk because the thermostat acts only at the 
boundary. In this subsection we derive an expression for the temperature of the reservoir similarly to how it has 
been done in Ref. [ p8| . 

If we assume equipartitioning of energy of particle and reservoir at the wall, we can define the reservoir temperature 
Tr indirectly via the velocity distribution of the particle at the moment of the collision denoted as Qmap- For sake of 
simplicity, here we do not explicitly consider the dependence of Tr on the position /3 of the colliding particle at the 
disk. An expression for the temperature of the reservoir can then be derived from the temperature in the bulk on the 
basis of the relation between the map density Qmap and the time-continuous density q in the bulk as given by Eq. (5) 
in Ref. |2|, 

g{v) ^ const.^^^^^ . (12) 

V 

The precise derivation of this equation can be found in Sect. IIIB2 of Ref. ||2^. To obtain the expressions for the velocity 
fluctuations parallel and perpendicular to the reservoir, the corresponding equations for the velocity distributions of 
the normal and tangential components u„ and Vt, resp ectively, have to be calculated. For this purpose, first the 
analogous equation for q{'j) corresponding to Eq. (|lj) must be derived. Knowing that in equilibrium £1(7) = 1 
because of symmetry and Qmapil) = cos 7 at the disk leads to 

e[^)^const.f^ (13) 
|cos(7)| 

Combining these two equations yields the full transformation 

COTlst . 

q{v)q{j) = — '--T7Qmap{v)gmaph) ■ 

i;|cos(7)| 

Changing to local Cartesian coordinates (vmVt) co-rotating with the position (3 at the disk and applying the trans- 
formation dvndvt = vdvd'y results in 

const 

Q{Vn)Q{vt) = '--T7Qmap{Vn)0mapivt) (14) 

ti|cos(7)| 

Noting that |u„| = v\ cos(7)| and matching the variables on both sides, Eq. (|lj) can be decomposed into 

Q{vt) = Qmap{vt) (15) 
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with —oc < vt < oo and 



const. 

ei^n) = -j T^QmapiVn) (16) 

\Vn\ 

with < ?;„ < oo. Before we come to the reservoir temperature definitions which are based on these densities, we 
remark that the disk which serves as the thermal reservoir is fixed and cannot recognize any current. In other words, 
only the kinetic energy of the particle in the fixed frame of the bulk 2Epf =< > is relevant for the interaction with 
the reservoir, and no average current needs to be subtracted. Defining now [. . .] as the average over Qmap Eq. ( [T5| ) 
implies for the tangential component < >= [v^] thus leading to the definition of Tj as 



Tt = 



V. 



4 

Analogously, the average over the map density corresponding to < t;^ > can be calculated from Eq. (Qq) to 



2 lZo'"n9{Vn)dv^ J^^KlQmapMdVn [\vn\ 



where the denominator is obtained from joint normalization over the ingoing and outgoing fiuxes. T„ is then defined 
as 

1 fl7,„l 

Tn — — 



with — oo < w„ < and < < oo. 

The total temperature of the reservoir T^. is consequently the average of Tj and T„, 

Tr = . (17) 

Tr can be calculated as an average over (3 or locally in a small interval A/3. We note that our result for T„ is slightly 
different to the result in Ref. |^8|, which strictly speaking is only valid if the in- and outgoing densities are symmetrical. 

This definition of the temperature is exact in equilibrium, however, in case of a nonequilibrium situation Eq. ( p^ ) 
and Eq. (|l3|) are not valid anymore. A more detailed analysis of these shortcomings leads to the conclusion that Tr 
calculated according to Eq. (|l^) will be greater than the real temperature of the reservoir for higher field strength 
prf . One would only obtain the real temperature of the reservoir if one would use the correct relation between map 
density and time-continuous density in nonequilibrium, and this is not known. In any case, a lower bound for the 
temperature of the reservoir which we denote as Tir can be calculated by only taking into account the velocity of the 
particle after a collision. 



V. EQUILIBRIUM 

The numerical calculation of the Lyapunov spectrum for the Lorcntz gas as thermostated by deterministic 
scattering according to the method presented in section || leads to the following result in equilibrium: {A} = 
{1.8695,0.0104,-0.0104,-1.8695}. These data appear to be at variance with the fact that in equilibrium two zero 
Lyapunov exponents have to exist, one associated with the direction of the flow, and a second one resulting from 
the conjugate pairing rule in equilibrium p^ . We have performed the following tests to detect the reason for this 
discrepancy: 

1. We have numerically calculated the Lyapunov exponents by investigating the dynamics of a trajectory and four 
satellite trajectories, i. e., for finite but small distances. The Lyapunov spectrum obtained by this method was 
the same. 

2. Changing parameters like the interdisk distance w, the parametric temperature T, the dimensionality of the 
reservoir 1^,0, the slope of the baker map, and replacing the baker map by more complicated two-dimensional 
maps like the cat map or the standard map pJ] did not improve the result. 
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3. We have followed the temporal evolution of two points on the same trajectory for about 20 collisions without 
orthonormalization and by choosing as initial conditions a) that the points are slightly displaced along the 
trajectory but have the same velocity, b) that the points have the same configuration space coordinates but 
slightly different velocities. We could then show that two neutral directions exist corresponding to a) the 
direction of the flow and b) to one direction perpendicular to the flow. 

The third test indicates that two zero Lyapunov exponents indeed exist. Thus, there must be a numerical problem 
because of standard Gram-Schmidt orthonormalization which establishes an orthonormal system of the tangent vectors 
on the basis of the most unstable direction. To cure that problem, we propose an alternative method to perform the 
periodic orthonormalization. This method establishes an orthonormal system of the tangent vectors starting from the 
existing neutral direction of the flow. Since we are introducing an additional constraint this way we call it constraint 
orthonormalization. It consists of the following steps: 

1 . Choose suitable initial conditions for the orthonormal system: The first tangent vector is situated in the direction 
of the flow, ^Fi = (vxo/voT^yQ/vo,0,0), and the other tangent vectors are orthonormal to it. 

2. At every orthonormalization the first tangent vector is forced to point in the direction of the flow, STi = 
ivx/v,Vy/v,0,0). This step corrects the very small deviations of the first tangent vector from the direction of 
the flow resulting from a collision with the disk, as will be explained in more detail below. 

3. The second, third and fourth tangent vector are orthonormalized again starting from the first one according to 
the method of Gram-Schmidt. 

The application of constraint orthonormalization leads to the following Lyapunov spectrum in equilibrium: {A} = 
{1.8695,0.0000,0.0000,-1.8695}, see also Fig. ^ Comparing these results to the previous ones obtained from the 
standard method shows that the Gram-Schmidt orthonormalization led to a wrong result only for the second and 
third Lyapunov exponent. The explanation for this numerical problem is as follows: In equilibrium the average energy 
transfer to the reservoir is zero. Still, at any collision energy is transfered either from the particle to the reservoir or 
in the opposite direction. According to Eq. (^) the phase space volume thus locally contracts or expands although the 
global phase space volume contraction is zero. However, as shown by Eq. (Q) the phase space contraction is intimately 
related to the corresponding (un)stable directions in phase space. Consequently, the local contractions and expansions 
at a collision change the orientation and the norm of the tangent vectors in a nontrivial way. The Gram-Schmidt 
orthonormalization reacts to these changes by turning the corresponding tangent vectors out of the previously neutral 
directions. The problem why the Gram-Schmidt procedure does not converge to the existing two neutral directions 
at least in the long time limit could not be completely resolved even by very detailed numerical investigations of 
the dynamics. Possibly some kind of resonance phenomenon between local phase space contraction and expansion 
at the collision and Gram-Schmidt orthonormalization after the collision leads to the corresponding tangent vectors 
adjusting themselves somewhat symmetrically around these two neutral directions p8| . 

In summary, constraint orthonormalization correctly yields a second zero Lyapunov exponent beside the zero 
Lyapunov exponent corresponding to the constrained direction of the flow. In agreement with the on average zero 
energy transfer between particle and reservoir the sum of the Lyapunov exponents and the global phase space volume 
contraction are zero. Furthermore, the Lyapunov exponents trivially fulfill the conjugate pairing rule related to the 
Hamiltonian character of the dynamics in equilibrium. The probability density in the Lorentz gas cell is uniform 
in equilibrium and, as a consequence, the Kaplan- Yorke dimension is equal to the dimension of the phase space, 
DjiY = d — 4:. In addition, the parametric temperature T, the temperature in the bulk Tf,, and the reservoir 
temperature Tr, are all equal, T = Th = Tr- 

We now turn to an even more detailed analysis of the dynamical instability of our model system by following ideas 
summarized in Refs. ^6| , ^ . If a dynamical system is ergodic, the Lyapunov exponents do not depend on the initial 
conditions of the tangent vectors, and thus they only yield information about the global instability. This implies that 
Eq. d^) provides no direct way to assess the local instability of the system at specific values of phase space variables 
like the angle of incidence at a disk 7 and the position of the colliding particle /3. In Refs. |3^,^, two slightly different 
ways have been proposed how to access information on local instabilities depending on these parameters. Here we 
use the approach proposed in Ref. which characterizes the local deformation of a typical tangent vector ST at the 
moment of a collision by introducing the quantity 

Ac(sin7,/3)=<lnffi>, (18) 

where the brackets indicate an average over all collisions in a respective small interval around (3 and/or 7. The 
physical motivation for defining this quantity is that any tangent vector quickly orients itself into the direction of 
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fastest growth. Accordingly, the full memory about the maximum instability of the system is contained in the 
orientation of the tangent vector thus representing a "needle" in phase space which very sensitively measures the local 
changes of the stability at a collision. Therefore, this quantity is a very sensitive function of 7 and /3. In Ref. ]39| ] 
this quantity has been called a local Lyapunov exponent, however, this term has also been used in the literature to 
indicate the dependence of the Lyapunov exponents Eq. (W)on initial conditions in case the dynamics is non-ergodic 
^2| . To avoid possible confusion, and by following Ref. [p6| where very related quantities have been defined, here 
we denote Ac as the local stretching rates of the system. Note that the clever and very simple definition by Eq. (fig) 
makes at least the maximum local stretching rate directly accessible to computer simulations. In contrast, in Ref. ||3^ ] 
the full spectrum of these rates has been defined in a proper co-moving coordinate system. This makes their definition 
more convenient in mathematical terms, but also less accessible for straightforward numerical computations. Both 
these different definitions are related via coordinate transformations |Q . Unfortunately, local stretching rates are not 
coordinate-invariant thus yielding different values depending on their precise definition, even in conjugate dynamical 
systems. 

Ac as a function of P for thcrmostating by deterministic scattering in comparison to elastic collisions is presented 
in Fig. |^(a) showing that the conventional hard disk Lorentz gas and our thermostated version of it share the 
same properties. The maxima/minima of Ac(/3) correspond to the directions of maximal/minimal distances between 
neighboring disks, respectively. Results for the conventional Lorentz gas in a more detailed view of phase space, i.e., 
Ac(sin(7)) for < 0.00001 as presented in Ref. have shown that the local stretching rate Ac is a singular function 
of sin(7) (see also Ref. |3^). Related results for Ac(sin(7)) for the Lorentz gas as thermostated by deterministic 
scattering are presented in Fig. ^. The curve in Fig. |^(b) looks qualitatively very similar to the curve in Fig. 1 
of 1^^. However, the numerical results for Ac(sin(7)) of our system are not sufficiently accurate |^ to study the 
existing discontinuities on a finer scale, as it has been done in Fig. 2 of Ref. [Q. To perform such investigations in 
a slightly more detailed way we looked at the refined, decomposed local stretching rate Acr^ =< ln(|(5r^|)/(|(5r2;|) > 
which characterizes the deformation of the x-component of a tangent vector only. The results for Acp^ are presented 
in Figs. H (d) and (e). Fig. ||(e) shows an enlarged sector of (d) where one can see a roughly symmetric profile 
composed of maxima and minima on a fine scale. The apparent symmetry of most of these peaks suggests that these 
oscillations are not due to numerical errors. We consider this as an indication that for the Lorentz gas as thermostated 
by deterministic scattering at least the refined local stretching rate Acr^ could be a singular function of sin(7). It 
may be somewhat surprising that such specific dynamical properties of the conventional, unthermostated Lorentz gas 
persist in our thermostated system as well. However, this leads to the conclusion that the geometric instability of 
the system is more important for these characteristics than the one resulting from the modifications related to our 
specific scattering mechanism. 



VI. NONEQUILIBRIUM 

In nonequilibrium we choose the electric field such that £x > 0, = 0. The field accelerates the particle, and energy 
is transfered on average to the disk resulting in a nonequilibrium steady state. As a consequence, the global phase 
space volume contraction given by Eq. (||) is negative. The detailed dependence of the Lyapunov spectrum on the 
field strength is shown in Fig. ^, where both results from the standard method as well as from the constraint method 
are presented. Only one zero Lyapunov exponent exists in nonequilibrium associated with the direction of the fiow. 
For higher field strength unconstraint Gram-Schmidt orthonormalization correctly turns the second tangent vector in 
the direction of the flow, A2 < 10^"* for Ex > 0.5. 

To obtain the correct Lyapunov exponents for Sx < 0.5 in nonequilibrium, we apply a suitably adjusted version 
of constraint orthonormalization as used in equilibrium: In order to achieve that two points on a trajectory stay on 
the same trajectory after a collision, their initial states and velocities have to be chosen such that the points have 
the same velocity at the moment of the collision. This condition leads to the components for the first tangent vector 
^Fi = {vxOTVyo,ex,0}- Note that ST is not normalized here. The other steps are then the same as in equilibrium. 

The Lyapunov spectrum as a function of the field strength obtained from constraint orthonormalization is also 
presented in Fig. ||. As in equilibrium, the results of the two methods differ only for the second and third Lyapunov 
exponent for Sx < 0.5. The differences for the third Lyapunov exponent are of the same size as the differences for 
the second Lyapunov exponent. The second Lyapunov exponent obtained by the new method is zero for all field 
strengths corresponding to the constrained tangent vector in the direction of the flow. The third Lyapunov exponent 
decreases with increasing field strength which is related to the dominant energy transfer in the direction from the 
particle to the disk. The dependence of the third and of the fourth Lyapunov exponent on the field strength appears 
to be a power law, which is a behavior that has also been observed for the Gaussian thermostated Lorentz gas for 
small enough field strength |Q. According to Pesins theorem [p2i, the only positive Lyapunov exponent is equal to 
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the Kolmogorov-Sinai entropy huSy 



T-KS 



Interestingly, its curve is nonmonotonous. For small field strength, the dynamics in configuration space appears to 
be dominated by the fact that the trajectory of the particle is getting adjusted in the direction of the field, and 
the Kolmogorov-Sinai entropy decreases. For higher field strength the increasingly disordered dynamics in velocity 
space related to an increase of the bulk temperature Tb seems to become more important, and the Kolmogorov-Sinai 
entropy increases. The same field dependence of the Kolmogorov-Sinai entropy has been observed in a shear fiow 
as thermostated by deterministic scattering [^o| . It would be interesting to know whether this is a general property 
appearing in field-driven system as thermostated by deterministic scattering. In contrast to this observation, the 
Kolmogorov-Sinai entropy monotonically decreases for the Gaussian thermostated Lorentz gas because the constraint 
of the bulk thermostat onto the dynamics increases with increasing field strength plp^ ]. Whether the irregularities 
on the fine scale in Fig. ^ are a property of the dynamics or whether they are numerical fiuctuations could not be 
decided on the basis of the present data. 

The sum of the Lyapunov exponents is negative and according to Eq. (^) equal to the phase space volume contraction 
P. As presented in Fig. ^(a), P decreases with increasing field strength. The density of the attractor remains phase 
space filling but shows a nonuniform and complicated structure as shown in the Poincare section in Fig. 4(a) of Ref. 
Therefore we can assume that the Hausdorff dimension Z?o is equal to the dimension of the phase space, Z?o = d = A, 
as is also the case for Gaussian thermostated periodic Lorentz gases. In contrast, the Kaplan- Yorke dimension Dky 
defined by Eq. (10) is not an integer anymore, as presented in Fig. ^c). This provides quantitative evidence for the 



fractal structure of the attractor according to the conjecture Dky = Di. 

Some of the conventional deterministic and time-reversible bulk thermostats fulfill the conjugate pairing rule saying 
that the Lyapunov exponents can be grouped into pairs such that A+ + A_ = const. [^,0. Fig. ^(b) shows that 
the conjugate pairing rule does not hold for our model. However, this does not come as a big surprise because it is 
well-known that even conventional thermostats do not exhibit conjugate pairing if thermostated at the boundaries 

The local expansion rate XdP) as defined in Eq. ( p^ ) is presented in Fig. ||(b) and can be compared with g{/3) shown 
in Fig. ||(a). One still recovers remnants of the periodic equihbrium distribution of Ac, see Fig. ^(a). However, they are 
strongly deformed by the anisotropy induced by the field, and the maxima and minima are much more pronounced. In 
contrast to equilibrium, there exist two absolute maxima, one around /? w 7r/6 and one around /? ~ 27r — 7r/6, and an 
absolute minimum around /? « tt. The maxima and minima of Ac(/3) occur just opposite to the maxima and minima 
of g{P). This is in agreement with the physical interpretation that a more unstable dynamics leads to a more dilute 
particle density in phase space. To extend the comparison, the temperature of the reservoir Tr as a function of /3 
calculated according to Eq. ( p!?] ) is presented in Fig. ||(c). The distribution of peaks in XdP) and Tr{P) is very similar. 
This might be related to the fact that both quantities illustrate somewhat irregular behavior: Ac characterizes the 
instability of the dynamics and Tr is equal to the mean kinetic energy of the degrees of freedom of the reservoir. The 
analogous /3-dependence of Ac and Tr points again to a close relation between dynamical system theory and statistical 
mechanics. At /3 « tt both the distributions of Ac and Tr show a more complicated structure. This is probably a 
consequence of the dynamics being directed parallel to the field resulting in the global minimum of Ac on a coarse 
scale, whereas for other /3 the dynamics is more chaotic. For more detailed views of the phase space in terms of the 
local stretching rate, in analogy to Fig. ^ in equilibrium, we could not get qualitative good results p| . Thus, whether 
Ac(/3) is a singular function in noncquilibriuni on a fine scale remains an open question. 

The dependence of Tr on the field strength according to the definition in Eq. (p^, in which Tr is averaged over /3, 
and the lower bound Tir as defined below this equation are shown in Fig. ^(d). In particular, the results for Tir{e) 
confirm that Tr is always greater in nonequilibrium than the parametric temperature T, Tr > T. 

More detailed information related to the deviations between reservoir temperature and parametric temperature are 
obtained by studying the map densities gmap{vt) and Qmap{vn) at a collision as represented in Fig. ^ The deviations 
between ingoing and outgoing densities in both cases are reminiscent of an average transfer of kinetic energy from 
particle to reservoir, as it is necessary to compensate the influx of energy caused by the electric field to obtain a 
nonequilibrium steady state. However, in case of the periodic Lorentz gas taking the thermodynamic limit leaves the 
system precisely as it is. Consequently, there is no thermodynamic way to get rid of the difference between ingoing 
and outgoing velocity distribution. But these differences are the dynamical reason why in nonequilibrium the reservoir 
temperature Tr is typically not equal to the parametric temperature T, because this would only be the case if both 
distributions would be converging to the (local) equilibrium distribution in the thermodynamic limit, as included in 
these figures. This aspect will become important for understanding our results on the relation between phase space 
contraction and entropy production below. 
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In Fig. 1^, the kinetic energy of the particle in the bulk averaged over /3, E^f =< > /2, is presented as a function 
of the distance d from the disk. The profile of Epf is inhomogeneous as expected. For c? — > 0, Epf should approach 
the temperature of the reservoir defined via equipartitioning of energy thus providing an alternative definition of the 
reservoir temperature based on the bulk dynamics in the inner ring around the disk. However, it is not possible to 
safely extrapolate to this limiting value on the basis of the present data Comparing Epf for rf ^ with Tr 

shows in particular that Tr > Epf for — excluding convergence via extrapolation, i.e., the assumptions made in 
the derivation of Eq. do not hold in nonequilibrium, but at least they yield a reasonable estimate. Note that 
the kinetic energy in this limit is still always greater than the lower bound for the reservoir temperature, Epf > T/^, 
which will be important for our following discussion of entropy production. 

The external driving force e performs work on the system and causes a macroscopic flow characterized by a positive 
conductivity tr > p7[ |. At the same time, work is transformed into heat and in turn removed by the thermostat 
leading to a positive thermodynamic entropy production. Starting from Eq. (pT|), the irreversible entropy production 
in the bulk dS is easily computed by defining the heat production dQ as the change of the kinetic energy of the 
particle in the bulk, dEpf/dt, and feeding in the bulk equations of motion Eq. This leads to the well-known 
expression of entropy production via Joule heating 

dS = f^L^^ . (19) 

1 r 

The numerical result for the field dependence of the thermodynamic entropy production according to this equation is 
presented in Fig. ||(a). 

On the other hand, as discussed above the heat produced in the bulk must leave as an outward flux across the 
walls absorbed by the thermal reservoir. Correspondingly, computing the average change of the kinetic energy during 
a free flight from the equations of motion yields 

< > - < t;'^ > 

ex < vx >= , (20) 

where the right hand side is just the average transfer of kinetic energy at a collision. Inserting this result into Eq. ( [l9| ) 
leads to 

dS = . (21) 

Comparing now Eq. (21) with Eq. (^ yields the important result that the identity between thermodynamic entropy 
production and phase space volume contraction does not hold for the Lorentz gas as therniostated by deterministic 
scattering. Instead, these two quantities just differ by the factor T/Tr{ex), 

To explicitly compare these two quantities the field dependence of —P is also presented in Fig. ^(a). As we have 
discussed above, there is some ambiguity in defining the reservoir temperature Tr, however, we emphasize that all our 
applied definitions and bounds lead to the result that —P and dS are inherently different in nonequilibrium. This is 
also clear from the fact how the thermostat works in our model, as explained above. 

As has been done in conventional thermostats, starting from Eq. (|2^) a relation between the electrical conductivity 
and the phase space volume contraction can now be established by using Eq. ( |l9| ) and replacing the average current 
according to the definition of the conductivity 



yielding 



TP -T^A 



This equation is formally identical to the Lyapunov sum rule obtained for the conventional thermostats. The only 
difference is the constant factor T which, for conventional thermostats, corresponds to the temperature of the reservoir. 
If the Lyapunov sum rule applies, it shows that macroscopic transport can be directly understood in terms of the 
microscopic dynamics characterized by the sum of the Lyapunov exponents. 
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However, we remark that the existence of a Lyapunov sum rule in thermostated systems being of the simple type 
as one above rather seems to be the exception than the rule: For example, a difference between phase space volume 
contraction and thermodynamic entropy production has also been obtained for a shear flow as thermostated by 
deterministic scattering |28[| . For this system the expressions for P and dS can be rather complicated. Consequently, 
the Lyapunov sum rule does not hold, and a similar relation has not been found in addition. Furthermore, already a 
variation of the Nose-Hoover and of the Gaussian thermostat did not lead to an identity between P and dS implying 
the invalidity of the Lyapunov sum rule as well, as discussed in Ref. pSp?] ]. 

In general, the relation between phase space volume contraction and thermodynamic entropy production, and 
the corresponding relation between transport coefRcient and Lyapunov exponents, will depend on the details of the 
microscopic energy transfer between particle and reservoir. Based on our studies in Refs. |4| p5|p7| , p8|| , we conclude 
that an identity between P and dS appears only to be valid for what might be called "ideal" thermostats meaning 
that energy is exchanged between subsystem and reservoir by sufficiently simple coupling rules as they are provided, 
for example by conventional Gaussian and Nose-Hoover thermostats. 



VII. CONCLUSIONS 



In this work we have numerically calculated the Lyapunov exponents for the Lorentz gas thermostated by deter- 
ministic scattering. The Gram-Schmidt orthonormalization, a fundamental ingredience of the standard method to 
calculate Lyapunov exponents, led to a wrong result for the Lyapunov spectrum by applying this thermostat. We 
modified this method by imposing an additional constraint, summarized as constraint orthonormalization, and found 
results which are in agreement with expectations from dynamical systems theory. We wish to remark that the phe- 
nomenon causing our numerical difficulties is reminiscent of an inelastic collision of a particle with a hard disk, as it 
is also modeled in granular materials by using restitution coefficients. Thus, applying constraint orthonormalization 
might be helpful for exactly computing Lyapunov spectra in low-dimensional systems of granular type as well. On the 
basis of the Lyapunov exponents further quantities have been calculated to characterize the nonequilibrium steady 
state. The comparison of the results obtained for thermostating by deterministic scattering with the ones known for 
conventional thermostats leads to the following conclusions: 

1. The sum of the Lyapunov exponents for thermostating by deterministic scattering is negative in nonequilibrium 
in agreement with the phase space volume contraction onto an attractor. For thermostating by deterministic scattering 
only one Lyapunov exponent is zero in nonequilibrium related to the direction of the flow. Similar results could be 
expected for the Nose-Hoover thermostated Lorentz gas where the calculation of the Lyapunov exponents have not 
yet been performed. In contrast, for the Gaussian thermostated Lorentz gas two Lyapunov exponents are zero in 
nonequilibrium because the thermostat keeps the kinetic energy of the particle strictly constant. 

2. The Kaplan- Yorke dimension calculated on the basis of the Lyapunov exponent is not an integer in nonequilibrium 
providing quantitative evidence that the attractor of thermostating by deterministic scattering in the periodic Lorentz 
gas exhibits a fractal structure analogous to the conventional bulk thermostats. 

3. The identity between thermodynamic entropy production and phase space volume contraction does not hold for 
thermostating by deterministic scattering. Instead, these two quantities differ by a field dependent factor. The reason 
for this difference is that the temperature of the reservoir of thermostating by deterministic scattering depends on 
the field strength, in contrast to Gaussian and Nose-Hoover thermostats. This result is important, since this identity 
was accepted up to now as a general characterization of nonequilibrium steady states generated by deterministic and 
time-reversible thermostats. 

4. Surprisingly, although there is no identity we could still establish a relation between conductivity and Lyapunov 
exponents for thermostating by deterministic scattering. This equation is formally identical to the Lyapunov sum 
rule for conventional thermostats. As far as we know, our model thus provides a first example of a system where 
there is no identity, but where nevertheless there is a simple relation between transport coefficients and dynamical 
instabilities similar to conventional thermostats. 

In summary, we find that the existence of fractal attractors in nonequilibrium steady states are common features 
which thermostating by deterministic scattering shares with conventional thermostats. Physically speaking, the fractal 
character reflects the extreme rarity of nonequilibrium states relative to equilibrium ones. To look for additional 
common properties of all deterministic and time-reversible thermostats remains an important question, which is 
intimately related to obtaining a general characterization of nonequilibrium steady states. Such a characterization 
might result in a more general relation between quantities of thermodynamic interest and the indicators of dynamical 
chaos at the microscopic level, from which the relations obtained for the thermostating mechanisms considered above 
could appear as special cases. 
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FIG. 1. (a) Elementary cell of the periodic Lorentz gas on a triangular lattice, (b) Definition of the relevant variables. 
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FIG. 2. Local stretching rate in equilibrium: (a) thermostating by deterministic scattering (solid curve), elastic collisions 
(dashed curve) for ah 7; (b)-(e) thermostating by deterministic scattering, (b) < 0.00001, (c) |/3 — 7r/2| < 0.00001; (d),(e) 
refined local stretching rate, (d) |/3 — 7r/2| < 0.00001, (e) enlarged sector of (d). 
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FIG. 4. Field dependence of various important quantities in noncquilibrium: (a) black: phase space contraction rate —P, 
white: irreversible entropy production dS; (b) conjugate pairs of Lyapunov exponents, black: Ai + A4, white: A2 + A3; (c) 
Kaplan- Yorke dimension Dky; (d) reservoir temperature Tr (black) and its lower bound Tir (white). 
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FIG. 5. Different quantities evaluated at the collision of a particle with the disk at field strength Sx = 0.5: (a) density at a 
collision q(I3); (b) local stretching rate Ac(/3), (c) reservoir temperature Tr{f3). 



16 



(a) 




(b) 1.0 



>^ 0.5 

Q. 



0.0 







1' 




V 




- /'/ 




/// 








/ 









1 2 



FIG. 6. Probability density at the moment of a collision for (a) tangential velocity vt and (b) normal component Vn at field 
strength = 0.5; solid curve: density before collision, dashed curve: density after collision, long dashed curve: density in 
equilibrium. 
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FIG. 7. Profile of the full kinetic energy Epf of the particle in the bulk, black Epf, white reservoir temperature Tr (a) field 
strength = 0.5, (b) = 1. 
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